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ABSTRACT 


In  airborne  anti-submarine  warfare  there  is  a  need  to 
more  accurately  determine  the  positions  of  sonobuoys  on  the 
surface  of  the  vater.   This  report  develops  two  algorithms 
which  employ  extended  Kalman  filters  to  determine  estimated 
position.   The  bearing  from  the  aircraft  to  the  sonobuoy  is 
the  primary  measurement.   Range  information  is  not  available. 
The  first  algorithm  is  a  six-state  filter  which  was  reduced 
from  the  13-state  system  developed  by  the  Orincon  Corporation 
Its  states  include  relative  jjosition,  relative  velocity,  and 
inertial  misalignments.   The  second  algorithm  includes  two 
cascaded  Kalman  filters.   The  primary  two-state  filter 
estimates  sonobuoy  position.   A  secondary  filter  estimates 
drift  from  information  obtained  from  the  primary  filter. 
Both  algorithms  successfully  estimated  sonobuoy  position  for 
simulated  aircraft  data.   The  effect  of  a ircraf t-to-sono bucy 
range,  the  frequency  of  measurement,  and  changes  in  altitude 
are  also  analyzed. 
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I.   INTRODUCTION 

The  P3  Orion  is  the  U.S.  Navy's  primary  long  range 
anti-submarine  warfare  (ASW)  aircraft.   It  is  outfitted  with 
equipment  which  allows  it  to  search  for,  locate,  and  track 
submarines.  The  aircraft  carries  a  Univac  digital  computer, 
the  CF  901,  v,hich  performs  much  of  the  navigational  and 
tactical  plotting  chores.   The  primary  sensor  used  by  this 
and  many  other  ASW  aircraft  is  an  airdropped  listening 
device  known  as  a  sonobuoy.   The  Orion  generally  deploys 
several  scnobuoys  (4  to  20)  in  patterns  which  can  cover  a 
thousand  square  miles  of  ocean  while  searching  for  the 
submarine.   Once  contact  has  been  made,  these  sonobuoys 
provide  information  which  locates  the  submarine.   The  target 
can  then  be  tracked  (or  attacked  if  required)  until  the 
Orion's  mission  is  complete. 

The  sonobuoy  is  dropped  from  the  aircraft  at  the 
geographical  location  designated  by  the  aircrev..   Once  in 
the  water  the  sonobuoy  floats  and  deploys  a  hydrophone  to 
depths  varying  from  50  to  500  feet.   The  information  picked 
up  is  transmitted  back  to  the  aircraft  where  it  is  analyzed. 
This  information  can  consist  of  the  relative  intensity  of 
target  noise,  bearings  and  sometimes  ranges  to  the  target, 
all  of  vhich  are  used  to  fix  the  current  position  of  the 
submarine.  In  order  to  maintain  close,  accurate  tracking  and 
be  able  to  launch  an  attack,  the  submarine's  position  must 
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te  accurate  to  within  several  hundred  yards.   The 
information  received  from  the  sonobuoys  has  some  error 
inherent  in  the  nature  of  the  measurements  made  ty  the 
hydrophone.   Sonobuoy  position  error  also  contributes 
significantly  to  the  submarine  tracking  inaccuracies. 

Historically,  the  positions  of  the  sonobuoys  were 
determined  and  updated  by  "mark-on-top"s .   This  required  the 
aircraft  to  home  on  the  transmitting  sonobuoy  until  the  buoy 
was  overflown.   At  that  instant  the  aircraft's  position  vas 
entered  into  the  on-board  computer  which  slewed  the  buoy  to 
this  updated  position.   After  many  of  these  updates  the 
computer  was  able  to  develop  a  bias  which  was  applied  to  the 
sonobuoy  positions  in  the  computer  effectively  allowing  them 
to  drift.  The  methoa  had  several  disadvantages.  Error  in  the 
updated  position  \*as  at  least  as  great  as  the  aircraft 
altitude  at  the  time  of  the  "mark"  which  could  vary  from  300 
to  20,000  feet.   The  accuracy  also  depended  en  the 
consistency  of  the  several  pilots  who  might  be  making  the 
mark-on-tcps  during  the  flight.  The  updating  was  not 
continuous  in  that  it  was  several  minutes  between 
consecutive  marks  on  one  buoy  at  best,  and  more  likely  30  to 
60  minutes.   Not  all  buoys  were  even  updated.   Furthermore, 
this  method  required  the  aircraft  to  overfly  the  submarine 
many  times  in  order  to  make  the  mark-on-tops.   This  should 
be  avoided. 

The  purpose  of  this  thesis  was  to  investigate  some 
alternative  methods  for  accurately  fixing  the  position  of 
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sonobuoys.  They  should  allow  the  aircraft  to  stand-off  from 
the  sonobucy  field  and  still  produce  more  accurate  fixing 
than  the  historical  method  provided.   The  Naval  Air 
Development  Center  (NADC)  at  Warminster,  PA.  had  already 
partially  developed  such  a  system.   This  thesis  was 
undertaken  in  support  of  their  work  but  v*as  conducted 
independently.  Their  system,  the  Sonobuoy  Reference  System 
(SRS),  was  already  installed  on  the  aircraft  and  had  the 
capability  of  measuring  the  relative  bearing  to  any 
transmitting  sonobuoy.   Additional  information  available  for 
use  included  aircraft  heading,  altitude,  and  airspeed  as 
well  as  doppler  velocity  and  drift  angle.   Also,  the 
aircraft's  Inertial  Navigational  System  (INS)  provided 
geographical  position  although  the  Schuler  cycle  and 
inertial  drifts  could  make  this  position  several  nautical 
miles  in  error.   On  the  other  hand,  aircraf t-to-scnobuoy 
range  and  sea  surface  drift  information  were  assumed  not  tc 
be  available.   An  attempt  was  made  to  determine  sonobuoy 
position  at  least  relative  to  the  aircraft  with  a  less 
accurate  geographical  position  as  a  secondary  objective. 
Salman  filtering  techniques  were  used  based  primarily  on 
measurements  of  bearing  from  the  aircraft  to  the  sonobucy. 
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II.   FILTEItS 

Kalman  filtering  is  a  recursive  technique  for  estimating 
the  state  of  a  system.   It  was  developed  in  the  1962's  by  R. 
E.  Kalman  and  improved  upon  previous  methods  by  Wiener  and 
others.   The  Wiener  filter  is  based  on  frequency  domain 
designs  which  are  statistically  optimal  but  are  only 
applicable  to  stationary  processes.   The  Kalman  filter  is 
based  on  state-space,  time  domain  formulations  and  is 
especially  suited  to  digital  computers.   From  a  simplistic, 
one-dimensional  point  of  view,  the  Kalman  filter  recursively 
averages  noisy  measurements  to  provide  a  more  precise 
estimate  of  the  actual  value. 

Assume  that  a  system  can  be  linearly  modeled  with  state 
equations  in  matrix  form  as 


X.  =  |,  X,   +  A.   U    +  T   V 
k    ~k-l  k4     k-i  k-i     k-t  k-1 


(1) 


where  X   represents  the  states  of  the  svstem  at  the  Kth 
k 

interval.  (Only  the  discrete  case  \*as  be  considered  in  this 

study.)  <5      is  the  transition  matrix  and  is  used  to 
"k 

"propagate"  the  system  from  X  tc  £+1.   U   represents  the 

k 

control  input  to  the  system  and  W   represents  white, 

k 

guassian  noise  with  0  mean  and  Q   variance,  written  N(3,Q  ) 

k  k 

Measurements  are  required  to  update  the  system  and  are 
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described  by  the  linear  matrix  equation 


Zk  -  \\   +   \ 


(2) 


where  Z,  is  the  measurement.   H,  describes  the  relationship 
k  k 

between  the  states  and  the  measurement  and  V,  is  measurement 

k 

noise  described  by  N(0,R  ). 

k 

The  discrete  Kalman  filter  equations  are: 
Propagate 


K  K  (+ 

k-i  k-i 


(3) 


?k(") 


X       p   (+)  (jjT   +  Q 
k-i  k-i      k-i     k- 


U) 


Update 


S     -  pk<->  < 


g  p  ( _ )  ff T  +  d 
nk  kv  ;  "k    nkj 


-1 


(5j 


Xk(+)  =  I  (-)  +  Gk 


Zk  "  W 


(6) 


P. (+)  ■ 
k 


I  -  G  H 
k  k 


V-' 


(?) 


vhere 


A 

X. 


=  estimate  of  the  state  X. 


covanance 
Salman  gain 
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A 

The  estimate  prior  to  the  measurement  is  denoted  by  XK(-) 
and  can  be  updated  after  the  measurement  to  a  new  estimate 

A 

denoted  XK(+).   This  notation  is  shown  in  figure  1.   The 
covariance  matrix  provides  a  statistical  measure  of  the 

A 

uncertainty  of  X.   Consider  a  2  x  2  covariance  matrix  where 
the  error  T  in  the  estimate  is  defined  as   x  -  x  : 


P   = 


E(x*  ) 


E(x.  lz  ) 


E(x,  x 


I  -""2 


E(x*  ) 


a,  a- 


(£) 


The  diagonal  elements  represent  the  mean  square  errors  of 
the  corresponding  state  variables  x   and  x2  .   The 
off-diagonal  elements  are  indicators  of  cross-correlation 
between  the  states.   The  Kalman  gain  J  is  an  optimal  gain 
choosen  so  as  to  minimize  the  sum  of  the  diagonal  terms  of 
the  covariance  matrix. 


Hk-i,Rk-i 


hl  a 


*k-iH 


pk-.'-i 


?k-.<*> 


Vi< 


£k«-l 


k,*k 

*k(*l 


*k-i,°k-i 


O 


*k.°k 


S> 


<k-i 


TIME  ► 


Eigure  1.   Salman  filter  notation 
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If  either  the  system  mcdel  or  the  measurement  model  is 
non-linear,  then  non-linear  filters  must  be  used.  Let  the 
non-linear  system  be  described  by 


\  -  ^K~r>  +  v,    ;  w,^n(0,q, ) 


(9. 


zk  =  h(v  +  vk 


I  Vk-K(0,Rk) 


(10) 


where  ©(X,  )  represents  the  non-linear  equations  cf  the 
r      k-i 

system.   They  can  be  linearized  cy   a  Taylor  series  expansion 
around  the  latest  estimate  of  the  state  tc  obtain  (j>  ( X ,)  . 
Likewise,  h(X  )  represents  the  non-linear  measurement 

A 

equations  and  by  a  Taylor  series  expanicns  yields  H(X.  ). 


(xk)  =  j*uk)  +  I(ik) 


X   -  X 
k     k 


11) 


h(X, 


h(Xk)  +  H(Xfc) 


Xk  "Xk 


(12) 


where 


fu )  =  — 


3_ 
dx 


J  x-x 


(13) 


H(X,  ) 

k 


a 

3x 


h(X, 


x--x 


14) 


Second  order  Taylor  series  terms  are  neglected.  Now,  the 
extended  Kalman  filter  can  be  implemented  with  the  following 
equat  ions : 
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Propagate 


X,  (-) 

k 


=  4  {* 


k-i 


(15) 


V"* 


-  iUvnt-))  ^.,(  +  )  i(x,  ,(- 


k-l 


•k-i 


(16) 


Update 


G   =  P  (-)  H(Xk( 


-)) 


H(X, .("))  P.  (-)  E(X.  (-)  )T- 


17! 


v+ 


v->  +  ^ 


zk  "  H(xfc(~) ) 


(18) 


¥ 


.(  +  )  = 


1  "  *kH(Ik(-))   Pk(- 


(19) 


Higher  order  filters  can  be  used  if  the  linearization 
errors  are  large.   The  second  order  Kalman  filter  employs 
one  more  term  in  the  Taylor  series  expansion  by  modifing  the 
update  equations  of  the  extended  Kalman  filter  to  account 
for  this  tern.   The  iterated  extended  Kalman  filter  uses  the 
same  equations  as  does  the  extended  Kalman  filter.   However, 
the  calculations  are  repeated,  each  time  linearizing  about 
the  most  recent  estimate,  until  there  is  little  further 
improvement  with  each  new  iteration.   The  iterated  extended 
filter  can  greatly  reduce  the  errors  due  to  non-linearities, 
mere  so  than  the  second  order  filter. 

Kalman  filters  should  be  based  en  correctly  modeled 
systems  and  accurate  noise  statistics  to  ensure  proper 
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performance.   This  is  not  always  possible  either  due  to 
ignorance  about  the  system  or  lack  cf  sufficient  statistical 
information.   A  filter  which  is  not  operating  properly  may 
diverge.   Apparent  divergence  describes  that  situation  where 
the  true  estimation  errors  are  larger  than  the  filter 
predicted  errors  although  they  are  bounded.   True  divergence 
is  characterized  by  errors  which  continue  to  grovi  with  time 
and  eventually  become  infinite.   These  divergence  phenomenon 
are  depicted  in  figure  2. 

There  are  several  ways  to  overcome  the  divergence 
problem  when  the  modeling  is  not  completely  accurate  such  as 
adding  fictitious  noise.   This  allows  the  filter  a  little 
more  freedom  to  adjust  to  whatever  modeling  inconsistencies 
may  exist,  but  makes  the  filter  estimate  appear  more 
erratic.   Another  method  *hich  helps  overcome  divergence  is 
finite  memory  filtering.   Since  Kalman  gains  tend  to  grow 
smaller  and  smaller  as  time  passes,  they  may  reach  a  point 


TRUE 


THEORETICAL 


TIME 
(a)  Apparent  Divergence 


TIME 
Ob)  True  Divergence 


Figure   2.    Divergence 
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where  new  measurement  information  has  no  effect  on  the 
estimates.   Finite  memory  filtering  effectively  eliminates 
old  data  which  is  no  longer  useful  by  keeping  the  gains 
significant.  This  is  sometimes  called  a  "moving  window  . 

Some  simplifing  techniques  used  in  Kalman  filters 
include  precomputed  gains.   Although  forfeiting  the  optimal 
Kalman  gains,  this  has  the  advantage  of  reducing 
computational  burden.   More  importantly,  the  gains  can  be 
controlled  to  overcome  modeling  weaknesses  if  necessary. 
Another  technique  which  can  be  used  when  more  than  one 
measurement  is  provided  to  a  filter  is  processing  them  one 
at  a  time.   This  avoids  taking  the  inverse  of  more  than  a 
scalar  when  computing  the  updated  covariance  matrix  P„(+). 
This  is  possible  if  the  simultaneous  measurements  are 
considered  to  be  taken  sequentially  over  a  zero  time  span. 
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III.   MODELS 

A.   GENERAL  APPROACHES  TO  SYSTEM  MODELING 

The  aircraf t-sonobuoy  system  must  be  modeled  in 
state-space  for  use  with  the  Kalman  filter.   There  are  at 
least  two  approaches  to  the  modeling  of  this  problem 
depending  on  the  point  of  view.  One  intuitive  approach  is  to 
assume  the  sonobuoy  drifts  at  a  constant  velocity  and  then 
use  aircraf t-to-sonobuoy  bearing  measurements  to  locate  the 
sonobuoy.   The  states  become  sonobuoy  position  and  sonobuoy 
velocity.   Unfortunately,  this  problem  is  not  observable. 
The  bearing  measurements  provide  only  information  about 
position?  there  is  no  rate  of  change  information  in  the 
bearings  themselves.   In  addition,  the  aircraft  must 
maintain  a  track  of  its  geographical  position  between 
updates  in  order  to  determine  the  next  expected  measurement 
for  the  Kalman  filter.   As  mentioned  before,  this  aircraft 
position  is  subject  to  non-linear  as  well  as  linear 
navigational  drifts  which  are  not  taken  into  account  in  this 
model.   However,  an  observable  system  can  be  obtained  by 
reducing  the  number  of  states  to  sonobuoy  position  only  and 
introducing  fictitious  process  noise  to  account  for  the 
drift.   This  noise  effectively  allows  sonobuoy  position  to 
update  so  as  to  keep  up  with  the  drift.   One  approach 
developed  by  this  thesis  is  a  variation  of  this  concept. 
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Another  approach  considers  only  the  relative  position  of 
the  sonotmoy  with  respect  to  the  aircraft.   Relative 
velocities  and  relative  accelerations  must  be  taken  into 
account  and  these  change  radically  as  the  aircraft  flies  in 
the  tactical  situation.   Sonobuoy  jjosition  is  not  obtained 
directly.   (When  only  the  word  "position"  is  used  it  will 
indicate  the  position  relative  to  an  earth  fixed  coordinate 
system,  such  as  latitude  and  longitude.   "Relative  position" 
will  always  mean  the  location  with  respect  to  an  aircraft 
fixed  coordinate  system.) 

There  are  some  other  considerations  which  should  be 
addressed.   The  sonobuoy  drift  is  generally  slow  (less  than 
5  NM/Hr  most  of  the  time)  and  constant.   It  is  not 
unreasonable  to  assume  that  the  entire  sonobuoy  field  drifts 
at  the  same  velocity.   Another  point  is  that  aircraft 
navigational  drift  can  not  be  distinguished  from  sonobuoy 
drift.  In  other  words,  the  drift  that  is  preceived  by  the 
aircraft  is  the  combination  of  sonobuoy  drift  and  aircraft 
navigational  drift.   If  this  navigational  drift  is  linear 
and  not  excessive  it  causes  few  problems.   dowever, 
non-linear  navigational  drifts,  such  as  the  Schuler  cycle, 
can  cause  large  errors. 
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B.   THE  COORDINATE  SYSTEM 

A  right-handed  coordinate  system  was  chosen  as  depicted 
in  figures  3  and  4.   This  is  a  departure  from  the  work  of 
Orincon  Corporation  which  is  described  in  the  next  section. 
This  system  allows  all  angles  to  be  measured  positive  in  the 
direction  they  are  normally  defined,  i.e.,  aircraft  heading 
measured  clockwise  from  north.   It  also  coincides  with  the 
usual  aerodynamic  coordinate  system. 


North 


East 


H   =  Aircraft  heading,  +  clockwise 


3 
Doi*n 

Eigure  3.  Earth  fixed 


Nose 


Aircraft  pitch,  +  up 
Aircraft  roll,  +  r.  wing  down 
Doppler  drift  angle,  +  r .  drift 
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C.   THE  SIX-STATE  SYSTEM 

Orincon  Corporation  of  La  Jolla,  CA.,  completed  a 
technical  report  in  December  of  1978  concerning  the  Sonobuoy 
Reference  System.   They  were  contracted  by  NADC  to 
investigate  the  effect  of  errors  in  the  aircraft 
navigational  system  on  the  performance  of  the  sonobuoy 
tracking  procedure.   As  a  result  of  their  v/ork  they 
developed  a  13-state  mathematical  model  which  successfully 
determined  sonobuoy  positions  in  spite  of  these  navigational 
errors.   The  complete  state  vector  was 


XT  = 


AX1  AX2  AV1  AV2  Xl  XZ    Vl  V2  *1  Xh    *3  %  m2 


(22) 


where 


Ax 
AV 
x 

V 

m 


relative  position  of  the  sonobuoy  from  the  aircraft; 

relative  velocity  cf  the  sonobuoy  from  the  aircraft; 

position  of  the  sonobuoy  on  the  ocean  surface; 

velocity  of  the  sonobuoy  on  the  ocean  surface; 

inertial  system's  platform  azimuth  alignment  errors 

inertial  system's  gyro  drift  rates  assumed  constant. 


(Subscripts  1,  2,  and  3  represent  North,  East,  and  dovn 
respectively.)  The  model  was  developed  using  state  equations 
which  described  sonobuoy  motion  relative  to  the  aircraft  and 
relative  tc  the  ocean's  surface.   These  were  derived  by 
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Orincon  frcm 


AX  =  X  -  (X1  +  E 
a 


(21 


and  other  mechanization  equations  to  obtain 


A*.  =  A? 


"1 
AX2 
Av 


1 


Av, 


1 

=  AV2 

2 
=  -R  j_  +  co  (  x  1 

■  ~R2  +  oo  U2 

=  v. 


-  i 


al  ~  A  V  "  ^3A2 


(Kg 


-  X 


a2 


-  A^2)  ~*3A!  "'^6 


Z2  =  V2 

7=0 
V„   =  0 


(22) 


where 


X1 


H  = 

g  = 

to  = 

A  = 


aircraft  inertial  navigation  system  position 

aircraft  inertial  navigation  system  error 

vector  from  the  center  of  the  earth  to  the  aircraft 

acceleration  due  to  gravity  (assumed  constant) 

Schuler  frequency  yg/R ' 

aircraft  accelerometer  outputs 


The  sonotuoys  vere  assumed  to  drift  at  a  constant  velocity 
and  the  aircraft  was  assumed  to  fly  at  a  constant  altitude, 
The  inertial  alignment  errors  were  described  by 


2k 


i>3   =     ij^n   cosx 


(23) 


inhere 

Q,   =   earth  rate 
X  =  latitude 

There  were  several  sources  of  measuremer. t  information. 
Aircraft  to  sonobuoy  bearing  v,as  available  from  the  aircraft 
interferometer  system  in  the  form  of  cosQ   vhere  8  is  the 
angle  between  the  direction  to  the  scncbuoy  and.  the  base 
line  of  the  antenna.   Doppler  velocity  and.  drift  angle  from 
the  aircraft  doppler  radar  system  provided  information  about 
Av   and  Av„.   The  only  other  relevant  source  of  measurement 
information  vould  have  come  from  JPS  (Global  Positioning 
System)  which  would  have  provided  error  free  aircraft 
position. 

Orincon  performed  a  numerical  observability  analysis  on 
this  system  based  on  the  observability  matrix  M,  given  by 


\-  I 


i=l 


£   H*  R_1H   $ 
k   k      k   k 
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If  this  matriz  is  positive  definite  for  H^2    then  the  system 
is  considered  observable.   In  addition,  an  eigenvalue 
analysis  provided  information  as  to  the  conditioning  cf  the 
system.   The  result  of  their  analysis  indicated  that  the 
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system  was  not  fully  observable  with  measurements  from  the 
interferometer  and  doppler  systems  alone.   Only  ten  states 
were  observable;  sonobuoy  position  most  likely  was  not 
observable  and  \p  ~,    m  ,,  m  2,  v.  and  v   were  weakly  observable 
at  test.   They  also  noted  that  "observability  is  reduced 
when  the  aircraft  pursues  a  straight  flight  path".   By 
removing  sonobuoy  positions  and  gyro  drifts  from  the  system 
they  found  the  observability  improved.   Further,  the  system 
became  fully  observable  if  sonobuoy  velocity  was  also 
removed  from  the  state  vector. 

The  Orincon  report  presented  the  results  of  a 
Monte-Carlo  simulation  on  the  13-state  system.   However, 
they  concluded  that  a  reduced  state  filter  made  up  of 


T 

X  - 


AX1  AX2  AV1  AV2  *1  *2 


(25) 


would  provide  a  workable  solution  to  the  problem  although 
they  noted  that  it  had  a  tendency  toward  divergence.   It  was 
partly  the  intention  of  this  report  to  further  the 
investigation  of  this  reduced  state  filter. 


26 


The  state  equations  for  the  reduced  six-state  filter  can 
be  obtained  directly  from  Equations  (22)  and  (23). 

Ax,  =  Av, 


Ax_  =  Av. 


Av      =  -R.  +    (j    {xl 

Av2    =  -I2  +   J(x2 

ty.    =  4>2  sinX 

i>2    =  -\pl  sinX 


xal    "   A21}    +   V 
Za2    "    Ax2>    "   *1* 


(25) 


where  ill  ,  m  ,  and  m   are  set  to  zero.   In  matrix  form 
3    1        2 


X  = 


[AV 

Ax„ 

AvL 

= 

AV2 

*1 

L  *J 

0        0 

1 

0 

0 

0 

0        0 

0 

1 

0 

0 

-co2    0 

0 

0 

0 

-g 

01             2 

0         _(JJ 

0 

0 

g 

u 

0 

0       0 

0 

0 

0       - 

-fisinX 

0        0 

0 

0 

f2sinX 

0 

Ax 


Ax. 


Av 


Av. 


0 

0 
n      2/       I 

-R  +  oil  -I  ,  , 

1  i   al 

-K-+  0)1  "I  -, 

2  2   a2 

0 
0 


which  can  be  abreviated 


X  = 


AA  i  k^ 


X  +  B 


(27) 


Expressed  in  discrete  time  the  state  equation  is 


\+l  =   »  Xk  +*k 


(26) 
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The  discrete  formulation  of  the  plant  matrix  can  be  obtained 
from 


i  -ti 


Si  -  A 


•H     At 
=  e 


(29) 


Developing  this 


si  -  A 


S        0 

1 

0  1 

1 

2 

0 

0        S 
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1  1 

1 
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0 

2     n 

OJ         0 

S 

0  ' 

1 

0 

g 

.          2 

0         CO 

0 

s  I 

"  _L 

~s 

0 

0   0   0   0'    s   ft  sin  A 

I 

0   0   0   0  '  -ftsinA    s 


(30) 


from  which  the  determinant  of  the  system  is 


si  -  A 


=  (s2  +  u2)2(s2  +  n2sin2\)2 


31) 


The  inverse  of  this  matrix  and  its  Laplace  transform  are 
presented  on  the  following  pages  as  Equations  (32)  and  (33) 
respectively. 
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The  discrete  formulation  of  the  continuous  matrix  S  is 


'k+l 


Ak  = 


eAt  dX  B 


where  A  =  At  -  t 


(34) 


Hence 


{-h-  +    (x.  -  x1.  )  }  (1  -  cosojAt 
to       1     a  1 

{-%■  +  (x,  -  x1,)  }  (1  -  coscoAt 
to       z     a/ 

{-—  +  (i»(x.  -  x1,  )  }  sincoAt 
1     al 


{- 


CO 


io(x„  -  x  „ ) }  sincoAt 
2     a2 


3 
3 


(35) 


The  noise  term  is  W  ~  N(0,Q).   The  variances  are  estimated 
to  te  .31  NM  for  relative  position,  1.0  (NM/Hr)  for  relative 
velocity,  and  .331  radians  for  inertial  misalignments.   As  a 
result,  the  following  constant  diagonal  matrix  is  used  for 
Q: 


Q  = 


.31 


.31  0 

1.0 

1.3 
0  .031 


.301 
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The  same  measurement  information  is  available  as  before 
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in  the  larger  system,  namely  interferometer  bearings  and 
doppler  velocity  and  drift  angle.   These  introduce  non- 
linearities  and  must  he  modeled  as  in  Equation  (10). 


\  =  h<V  +  vk 


A  Taylor  series  expansion  as  in  Equation  (11)  about  the 
latest  estimate  of  sonobuoy  position  provides  a  method  of 
linearization.   As  a  result 


K  = 


3  Ax. 


3  z 


3  Ax. 


3  z 


3  Av. 


*■-* 


3  z 


3  Av, 


3  z 


3^ 


3  z 


3  i>, 


(3?) 


-i 


To  determine  z  and  its  partial  derivatives,  recall  that 


z   =  cosi 


=  AX  •  R# 


(35 


where  AX  is  the  unit  vector  representation  of  a  line  from 
the  aircraft  to  the  sonobuoy  and  RB  is  the  unit  vector 
describing  the  base  line  between  the  antenna  pairs  of  the 
in terf erome ter. 


z  .  =  cos  9 


v--x 


Ax. 


Ax 


RB 


Ax, 


Ax 


RB 


Ax. 

— -  RB„ 

Ax 


Ax  RB.  +  Ax„RB   +  h  RB 
(AxL2  +  Ax22  +  \})h 


where  h  is  the  aircraft  altitude 


(39) 
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inhere 


"  9zl 

3zi 

La  Ax. 

3  Ax2 

2      0      0 


(40) 


A-t 


a«x 


a  a. 


3  Ax, 


•(z 
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X  =  X 


(Ax.+Ax.+h    ) 
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(Ax.+Ax    +h    ) 


HBr 


R32~ 


Ax/AXjRB,  +A^IB2  +HRBP 
(AxJ    +    AxN^    +    h2) 

Axz(AxRB;  +  AxRBz  +hRBp 
(AXj    +    Ax2    +    h2) 


(41) 


(42) 


The  doppler  velocity,  V,,  is  a  measurement  of  aircraft 

d 

speed  relative  to  the  ocean  surface  and  doppler  drift  angle, 
3,  is  a  measure  of  the  angle  between  aircraft  track  and 
aircraft  heading.   These  measurements  can  he  related  to  the 
state  vector  by  assuming  the  sonobuoys  drift  at  the  same 
velocity  as  the  ocean  surface.   If  this  is  true  then 


z2  -  Vd  =  (Av2  +  Av2)% 


(43) 


Expanding  in  a  Taylor  series  yields 


z2  =  (AvA2  +  AS2)1* 


(44 
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a  Av. 
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(47: 
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Similarly, 


z  -  B  =  Track  -  Heading 
3         -i/Av„\ 


=  tan 


-  H 


vAv 


or   since    the   aircraft   heading,   H    ,    is   deterministic 


(46) 
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\Av 


1' 


2      0 


a  z. 


3  Av. 


3  z 


3  Av. 


x=* 


A^X 


where 

3z3 

3  Av 

3z3 

3  Av2 

Finally, 

-Av, 


jr-x 


AVj    +    Av2 


A 

Av. 


Av      +    Av 


x=x 


0       2 


(49; 


(50) 


(51) 


(52) 


Z    = 


CCS0 
V 


H    -    8 

e 


H    = 


3  (cos  6)      3  (cos  e) 
x,  x„ 


0 


0 
2 


0  0 

3(vrf)  3  (v.j  ) 

3Av!  3Av2 

9  0fe+8)  3(Htf+8) 


3Av, 


3Av„ 


0  0 

0  0 

0  0 


(54) 


3^ 


R   =  .0007  rad 
R  =  1.0  (NM/Hr) 

R3  =  .01  deg 


The  measurement  noise,  V,,  is  described  by  N(0,R).   It 
is  assumed  that  the  measurement  noise  R  and  the  plant  noise 
Q  are  white,  gaussian,  and  uncorrelated .   Based  on 
information  obtained  from  NACC,  R  was  assigned  the  following 

values  : 

( interferometer) 

(doppler  velocity)  (55) 

(floppier  drift  angle) 
It  should  be  pointed  out  that  the  doppler  measurements 
are  not  always  available.   The  doppler  system  freezes  the 
last  value  of  groundspeed  and  drift  angle  anytime  the 
doppler  radar  is  not  receiving  good  information.   This 
happens  whenever  the  aircraft  is  above  50<30'  altitude,  or 
can  happen  when  flying  above  a  cloud  cover  or  ^hen  the  sea 
surface  is  too  smooth.   As  a  result,  the  doppler  could  be 
inaccurate  during  a  significant  portion  of  the  flight. 
The  six-state  system  was  programmed  and  tested  for 
several  conditions.   Flow  charts  for  significant  portions  of 
the  program  and  the  entire  program  listing  are  presented  in 
Appendix  A.   The  main  program  has  two  purposes:   first,  to 
control  the  simulation;  and  second,  to  read  and  prepare 
measurement  data  for  the  subroutine  FILTER.   This  subroutine 
performs  all  the  Salman  filter  computations  as  specified  in 
the  above  equations.   The  Carlson  square  root  technique  is 
used  to  ensure  a  positive  definite  covariance  matrix.   Also, 
in  order  to  avoid  inverting  a  matriz  in  the  calculations  for 
the  Kalman  gain,  the  three  measurements  are  processed  one  at 
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a  time.   The  results  of  the  simulation  are  presented  in 
section  IV. 3. 
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D.   THE  TWO-STATE  SYSTEM 

An  attempt  was  made  earlier  on  in  this  research  tc  model 
sonobuoy  motion  as 


X,   =  X   +  U.  At 
k-H     k     k 


U,   =  U, 
k+l     k 


(56 


with  a  state  vector  made  up  of  x    ,  x  ,  u .,  and  u  .   This 
system  proved  to  be  unobservable  although  the  simple 
approach  was  appealing.   Since  sonobuoy  drift  rates  are 
generally  slow,  another  attempt  was  made  to  model  the  system 
without  velocity.   The  equations  reduced  to 


X,   =  X,  +  w, 

k+l     k     k 


(57! 


With  only  x   and  z   as  states,  this  system  was  observable. 
3y  describing  large  statistical  values  of  system  noise  there 
i*ould  hopefully  be  enough  freedom  in  the  "update"  tc 
compensate  for  any  sonobuoy  drifts.   Unfortunately,  this 
approach  was  not  completely  successful!  hovever,  it  vas 
observed  that  the  general  direction  of  sonobuoy  drift  kas 
correct.   As  a  result  of  these  investigations  a  technique  of 
cascading  Kalman  filters  was  used  to  solve  the  tracking 
problem.   The  intent  was  to  use  results  from  one  filter  as 
measurement  for  a  second  filter.   Then,  the  result  cf  the 
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second  filter  could  be  used  as  a  deterministic  control  input 
to  better  propagate  the  first  filter. 

Ti*o  Kalman  filters  are  used  in  this  approach.  Figure  5 
shows  how  they  are  related  (or  cascaded).  The  first  filter 
models  soncbuoy  position  and  is  based  on  the  discrete  state 
equation 


X,   =  X  +  U,  At  +  W, 

k      k     k       k 


(58) 


The  state  vector  consists  of  the  sonobuoy  positions,  Xj  and 
X2 •   U   is  a  deterministic  input  accounting  for  soncbuoy 
drift  and  W   is  system  noise.   The  interferometer 

fcC 

measurement,  cosQ,  is  used  to  update  this  filter. 
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Figure  5.  Cascaded  Kalman  filters 
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The  second  filter  models  the  sonobuoy  drift  with  the  state 
equat  ion 


k+l     k 


(59) 


As  the  sonobuoy's  position  changes  with  each  update  of  the 
position  filter,  the  velocity  required  to  move  from  X   to  X 
in  the  time  interval  At  can  be  calculated.   This  value  is 
used  as  a  measurement  for  the  drift  filter.   All  sonotuoys 
contribute  to  this  filter  which  produces  one  overall 
estimate  of  drift  for  the  sonobuoy  field.   This  value  is 
then  used  as  an  input  to  the  position  filter  when  the  next 
interferometer  measurement  occurs. 

The  technique  is  understandably  sensitive  tc  the  output 
of  the  drift  filter.   If  the  estimated  drift  is  too  much  in 
error,  the  estimated  positions  will  be  affected.   Divergence 
could  result.   7or  this  reason,  a  simplified  Salman  filter 
is  used  to  estimate  the  drift.   It  is  of  the  form 


U,   »  U  + 


\{\  '  \ 


ed  > 


The  gain  is  precomputed  from  l/(n+l)  where  n  is  the  total 
number  of  measurements  and  is  limited  to  a  value  at  least  as 
great  as  £t/(S4  min.).   This  limit  is  reached  after  the 
first  S4  minutes  of  tracKin-  time.   From  that  point  en  the 
gain  remains  constant  effectively  averaging  the  drift  over 
one  period  of  the  Schuler  cycle.   It  must  be  understood  that 
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the  drift  this  filter  estimates  consists  net  only  of 
sonobuoy  drift  but  aircraft  drift  as  well  and  is,  in  fact, 
the  vector  sum  of  these.   An  estimate  which  attempts  to 
follow  the  changing  velocities  through  the  Schuler  cycle 
tends  to  cause  the  position  filters  to  diverge.   This  might 
be  the  case  if  larger  values  of  gain  are  used,  averaging 
only  the  last  few  measurements  of  drift.   Therefore,  many 
measurements  of  drift  are  averaged  which  provides  an 
estimate  for  linear  drift  only.   The  lower  limit  of  gain 
prevents  new  drift  measurements  from  being  ignored  and  also 
provides  for  flexibility  in  estimating  the  drift. 

Unlike  the  six-state  system,  the  interferometer 
measurement  is  the  only  measurement  required  by  the  position 
filter.   The  measurement  equation  is  essentially  the  same  as 
in  the  preceeding  section  where  z  is  determined,  as  shown  in 
Equation  (39).   Since  there  is  only  one  term  in  z,  the  H 
matrix  is  a  row  vector  consisting  of 


n    = 


3  z 

3x 


3  z 


3  x. 


(61) 


where  the  terms  are  given  in  Equations  (41)  and  (42).   It  is 
necessary  to  determine  the  ai rcraf t-to-sonobuoy  range  for 
these  values  from 


AX  =  X  -  X 


(62, 


A 

where  X  is  the  estimated  soncbuoy  position.   X1  represents 
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the  aircraft  position  and  is  considered  deterministic. 
Undoubtedly,  X   is  in  error  due  to  aircraft  navigational 
drifts  and  therefore  causes  errors  in  the  estimated 
geographical  position  of  the  sonobuoy.   But,  since  the 
measurement  of  bearing  is  based  on  relative  positions,  the 
estimate  is  relatively  correct. 

An  observability  analysis  on  the  2-state  position  filter 
was  performed  again  using 


N 


■  li 

k=l 


k  Ek  *'\    ^k 


(63) 


Subtitution  into  this  formula  with  N=2  yields 


M      = 


f3Z 
13  x 


Idz 
l3x, 


/3  2 

3x, 


f  3  Z      I  Id  z 


\dx,  /.\3x 


Idz 


\dx. 


v3x. 


t+i 


'3z_ 

<3x, 


/3z   \  /a  z 


V3xi"i^ 


3x 


21+) 


i+A 


3z 

,3x . 


I+l 


'3  2 


13  x 


2A+I 


(64) 


Expanding , 


M      = 


'32 


\3x 


lA+l 


/3z 
Ux. 


'3z 


(3z 


'3  2 


^xU+iV9x2A+iV3vxiA 


'3  2 

i3x. 


'3  2 

\3x, 


'd  Z 


\3x 


2/ 


32 


f3z 


V3xli+lV3x2i  +  l 


(65) 


*H 


This  matrix  is  positive  definite  for  all  except  the  case 
when 


3x 


1/ 


3z 


3x 


and 


14+1 


3z 

3x. 


3z 
3x, 


21+1 


This  occurs  whenever  the  relative  bearing  from  the  aircraft 
to  the  sonobuoy  is  not  changing  (i.e.  when  flying  directly 
toward  or  away  from  the  sonobuoy).   This  is  intuitively 
correct  since  two  or  more  bearings  must  cross  in  order  to 
determine  a  position. 

The  two-state  system  was  programmed  and  tested  to 
determine  the  usefulness  of  this  simplified  filter.   The 
majority  of  the  program  is  the  same  as  the  one  used  for  the 
six-state  system.   The  significant  differences  occur  in 
subroutine  FILTER.   The  Carlson  square  root  technique  is 
used  again  but  only  one  measurement  instead  of  three  is 
processed.   The  drift  filter  is  programmed  in  this 
subroutine  along  with  the  position  filter.   A  flow  chart  and 
program  listing  are  ^resected  in  Appendix  A. 
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IV.  ANAYLSIb 

A.   THE  SIMULATION 

Actual  data  from  the  Orion  was  not  easily  obtainable  for 
this  research.   Consequently,  a  computer  program  had  to  be 
written  to  generate  the  information  required  by  these 
algorithms.   A  flow  chart  in  Appendix  3  describes  the 
program.   An  aircraft  track  was  created  by  alternating  lines 
and  curves  of  various  lengths  and  then  a  determination  of 
noise  free  measurements  was  made  as  the  track  was  flown. 
(Measurement  noise  was  added  later  during  the  simulations.) 
Sonotuoys  i*ere  allowed  to  drift  at  a  constant  velocity,  and 
the  aircraft's  navigational  drift  was  modeled  with  a 
constant  velocity  and  Schuler  cycle  variations  as 


X1  =  X   +  K  At  +  A  sinoot 
a     a     x       x 


(66) 


inhere 


X   =  aircraft  inertial  position 


X   =  aircraft  true  position 
a  r 


ui  =  Schuler  frequency  -^/g/R ' 
The  constant  drift  rate,  K  ,  was  found  by  NAOC  to  have  d 
mean  and  a  standard  deviation  of  2.5  NM/HH.   Likewise,  the 
amplitude  of  the  Schuler  cycle,  A  ,  had  2    mean  and  a 
standard  deviation  of  .5  NM.   The  effect  of  wind  was  also 
available  in  the  program  but  was  never  included  for  analysis 
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It  v.as  obvious  that  both  algorithms  vere  extremely 
sensitive  to  the  aircraft  flight  path.   Since  there  is  no 
"typical"  flight  path  for  an  Orion  during  its  on-station 
period,  several  uniform  patterns  were  selected  which  would 
provide  meaningful  information  on  eacn  algorithm's 
performance.   These  base  parameters  were  chosen: 

airspeed     180  kts 

altitude     3000  ft 

range        15  NM  ai  rcraf  t-to-s  onobuoy 

frequency    20  sec  between  measurements  on  one  sonobuoy 
A  circular  pattern  was  flown  around  a  sonobuoy  at  a  range  of 
about  15  hiv.   Initial  sonobuoy  placement  was  at  (1,1)  with 
the  aircraft  flying  clockwise  starting  at  (15,0).   This  is 
shown  in  Figure  6  which  depicts  a  portion  of  the  aircraft's 
true  track.   Sonobuoy  drift  was  west  at  1  N'M/HRj  aircraft 
navigational  drift  was  south  at  1  NM/HB  and  included  a 
Schuler  cycle  (amplitude  .5  NM)  when  specified.   (This  drift 
would  be  indicated  in  Figures  12  thru  30  by  the  following 
notation:  DR=(a/c  150-1  +  Schuler,  b  270-1)  .)   The  sonobuoy 
remained  near  the  center  of  the  pattern  allowing  the  range 
to  remain  relatively  constant.   This  was  desirable  since 
this  basic  flight  path  was  used  to  test  the  response  of  the 
algorithms  to  variations  in  altitude  and  rar.-e.   Similarly, 
a  square  pattern  >as  flown  counter-clockwise  around  a 


(x,  ,Xj  )  measured  in  nautical  miles  from  an  arbitrary 
origin  near  the  aircraft's  starting  point. 
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sonotuoy  at  a  range  of  15  MM.   Ir.  this  case  the  aircraft 
began  on-top  the  initial  position  of  the  sonobuoy  at  (-3,-3) 
and  preceeded  to  fly  the  track  shov.n  in  Figure  7.   The 
sonotuoy  was  allowed  to  drift  to  the  east  at  2  NM/HR  and  the 
aircraft's  navigational  system  to  the  north  at  2.5  NM/HR. 
This  pattern  was  used  eventually  to  analyze  performance  at 
different  measuring  frequencies.   Figure  8  shows  the  output 
of  the  aircraft's  navigational  system  when  it  has  a 
northerly  drift  of  2.5  NM/HR  while  the  aircraft  is  flying 
the  square  pattern.   Figure  9  shows  this  same  pattern  when  a 
Schuler  cycle  with  an  amplitude  of  .5  NM  is  also  present. 
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Figure  6.   Aircraft  track  for  the  15  NM  circular  pattern 

with  initial  sonobuoy  location  and  direction  of  drift  indicated 
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Figure  7.   Aircraft  track  for  the  IS  NM  square  pattern 

with  initial  sonobuoy  location  and  direction  of  drift  indicated. 
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Figure  8.   Aircraft's  navigational  output  for  the  15  NM 
square;  navigational  drift  is  2.5  NM/Hr  to  the  north. 
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Figure  9.   Aircraft's  navigational  output  for  the  15  NM  square; 
navigational  drift  is  2.5  NM/Hr  north  with  a  Schuler  cycle. 


^9 


These  patterns  were  used  to  generate  the  data  vhich 
included  sonotuoy  number,  time  of  measurement  and 
measurement,  aircraft  latitude  and  longitude,  altitude, 
heading,  pitch  angle,  roll  angle,  N-S  acceleration,  E-W 
acceleration,  and  the  antenna  used  for  the  measurement. 
Also  determined,  tut  used  only  for  error  analysis,  were 
sonobuoy  position  and  true  aircraft  position.   This 
information  was  computed  for  the  entire  simulation  period 
and  stored  in  a  file  to  be  used  when  needed. 

Initial  estimates  of  the  state  and  the  square  root  of 
the  covariance  matrices  were  obtained  as  follows: 
Six-state 


X(0)  = 


x,(i3)-i  (0) 


zi(0)-2ai(0)  +  .5 

-1L  cosK 
-V,  sinH 

0 

0 


S  (0) 
D 


0  . 


.2025 
.0025 


(6? 


Two-state 


X(0)  = 


*l  (0) 
x2  ( 0 )  + .  5 


5D(0)  = 


r 


Alt 
Alt 


(68) 


where  S  (0)  indicates  the  diagonal  elements  of  5(0 
D 
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For  the  circular  pattern 


xl(d)    =   1. 
x2(0)  =  1. 


x  (0)  =  15. 
al 

X  (0)  =0. 
a2 


and  for  the  square  pattern 


xL(0)  =  -3. 
x2(0)  =  -3. 


x  (0)  =  -3, 
al 

x  (0)  =  -3. 
a2 


The  initial  values  of  the  covariance  were  chosen  so  as  to 
describe  the  errors  in  the  state  associated  with  dropping  a 
sonobuoy  from  an  aircraft.   The  two-state  system  operated 
tetter  vith  a  lower  initial  variance  than  did  the  six-state 
system . 

Simulations  were  run  for  a  nominal  period  of  two  hours. 
Each  unit  of  data  was  sequentially  read  into  the  simulation 
program  from  the  storage  file.  Measurement  noise  was  added 
with  the  following  normal  distributions: 

ccsQ  N(0,.0007)      in  radians 

doppler  velocity      N( 0,1.0)       in  NM/EH 
doppler  drift  angle   N(0,.01)       in  degrees 
These  values  agree  with  the  ones  chosen  for  H  as  described 
previously.   A  random  number  generator  using  an  initial  seed 
was  used  to  create  the  noise  from  the  proper  distributions. 
After  the  estimate  was  made,  errors  were  measured  ard 
manipulated  in  subroutine  RESULT.   This  cycle  was  repeated 
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until  the  simulation  was  complete.   (Refer  again  to  Appendix 
A  for  a  flow  diagram  of  the  simulation  process.) 

Information  vas  gathered  and  is  presented  primarily  in 
two  ways.   First,  for  both  the  square  and  circular  patterns 
the  results  of  the  first  run  of  the  simulation  are  plotted. 
(Figure  13  is  an  example.)   The  top  portion  of  the  results 
shows  the  estimated  positions  of  the  sonobuoy  on  the 
aircraft's  navigational  plot.   It  should  be  noted  that  these 
are  not  true  positions  since  inertial  errors  may  be  causing 
this  plot  to  drift.   The  bottom  portion  shows  the 
North-South  and  East-West  errors  in  tne  estimated  positions. 
These  errors  are  the  differences  between  the  true  and  the 
estimated  relative  positions  of  the  sonobuoy.   They  are 
depicted  in  Figure  10  and  were  computed  from 


AX,  =  AX,  -  AX, 
k      k      k 

AX,  =  (X  -  X1), 
k  a  k 


(X  -  X  ), 
a  k 


(69) 
(70) 


where  AX,  is  the  relative  error  at  time  k.   The  true  aircraft 
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Figure  10.   The  effect  of  the  aircraft  navigational  plot 
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position  is  denoted  by  X   and  the  true  geographical  position 
of  the  sonobuoy  is  denoted  by  X.   These  values  are  known 
from  the  simulated  data.   X  is  the  filter's  estimate  of 
position  and  X   is  a  deterministic  input  from  the  aircraft's 

cL 

navigational  system.   The  errors  indicated  in  these  plots 

are  different  from  those  perceived  in  the  upper  plot 

whenever  the  aircraft's  navigational  picture  is  drifting. 

The  positive  and  negative  values  of  the  square  root  of  the 

covariance  are  also  shown  on  these  plots  as  solid  lines. 

Specifically, 

N-S  error  =  ±-/P  ' 

E-W  error  =  ±  -/? 

Y      22 

(Only  jjosition  errors  are  analyzed  in  this  report.) 

The  second  way  in  which  the  information  was  gathered 
consisted  of  a  shortened  Monte-Carlo  simulation.   7or  each 
scenerio  under  study,  20  two-hour  simulations  were  run,  each 
with  new  values  of  measurement  noise  provided  by  the  random 
number  generator  from  their  respective  distributions.   The 
objective  was  to  compare  the  RMS  errors  predicted  by  the 
covariance  matrix  in  the  filters  to  the  actual  RMS  errors 
observed  in  the  simulations.   Three  RMS  statistics  v.ere 
collected  and  plotted  versus  time  as  follows: 

1.  The  square  root  of  the  covariance  computed  by  the 

Kalman  filter  is  represented  by  a  solid  line  on  the 

■clots.   It  was  commuted  from 


-LVp   +lYp   ,  n=20, 
n^  11   n  <— '  22 


k  =  1,2, 


(71 
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This  standard  deviation  is  the  filter's  estimate  of 
its  accuracy. 
2.  The  mean  and  standard  deviation  of  the  error  in  the 
estimated  sonobuoy  position  relative  to  the  aircraft 
was  determined  by 


y    =   -XA^ki'    n  =  20,  k  =  1,2,3...        (72) 
ni=l 


a2   =   -V^AX2  -  p.2,   n  =  20,  k  =  1,2,3. .. 
x.     nZ_i   ki    k 
k     1=1 


This  standard  deviation  is  represented  by  an  'X'  on 
the  plot  and  is  a  measure  of  the  variability  in  the 
filter's  estimated  .position. 
3.  The  RMS  value  of  the  actual  relative  error  was  also 
also  determined.   It  was  computed  as 


Li 

n  -^ — i 


i=l 


ki 


n  =  20,  k  =  1,2,3. 


(74 


and  is  represented  by  an  '0'  on  the  plot.   This  value 
is  a  measure  of  the  error  which  occurred  between  the 
actual  relative  position  of  the  sonobuoy  and  the 
estimated  one.   (See  Figure  11.) 
These  values  are  plotted  over  the  two  hour  simulation  period 
(k  =  1,2,3...)  and  provide  a  measure  of  the  accuracy 
obtained  by  each  of  the  algorithms. 


5k 


-S 36- 


Gx— *~ 
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B.   THE  RESULTS 

Orincon  Corporation  performed  much  of  their  analysis 
using  the  small  circular  pattern  shovn  in  Figure  12.   The 
aircraft  sped  around  this  track  at  36£  knots  taking  a 
measurement  every  20  seconds.   Altitude  was  constant  at  3000 
feet.   The  initial  position  cf  the  sonobuoy  was  at  (2,2)  and 
it  drifted  on  a  heading  of  345  degrees  at  5.5  N'M/HH  while 
the  aircraft's  navigational  system  drifted  in  the  opposite 
direction  at  5.2  NM/HR.   A  Schuler  cycle  was  alsc 
superimposed  on  this  drift. 

The  initial  estimates  of  sonotuoy  position  and  variances 
used  by  the  algorithms  were  different  for  this  particular 
pattern  ir.  order  to  coincide  with  the  Orincon  simulation  as 
much  as  possible.   They  were 

Two-s  tate 


M)   = 


0 


sD(0)  = 


2.] 


(74. 


Six-state 

-5.4 

2 . 

0 

2  . 

X<3)  - 

-Va  cosH 
-VasinH 

SD(0)  - 

5. 

0 

.  0025 

0 

.0025 

(75) 
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Vfi   =   2.8   M/S 


Figure    12.      Orincon's    simulation   flight   path. 
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Initially,  both  the  six-state  and  the  two-state 
algorithms  were  tested  using  this  flight  path  and 
comparisons  >ere  made  with  the  results  cf  Orincon's  13-state 
system.   It  was  discovered  that  the  six-state  system's 
results  did  not  change  as  a  consequence  of  vhether  or  not  a 
Schuler  cycle  existed  in  the  navigational  system. 
Therefore,  the  six-state  system  was  run  only  with  the 
Schuler  cycle  active.   The  two-state  system  was  run  with  and 
without  a  Schuler  cycle. 

The  results  of  the  first  run  of  the  simulation  are  shown 
in  Figures  13,  14,  and  15.   The  sawtooth  shape  of  the 
covariance  (solid  lines)  is  the  result  of  the  circular 
pattern  and  is  caused  in  two  ways.   Primarily,  since  the 
soncbuoy  remains  closer  to  one  side  of  the  circle,  there  is 
a  minimum  value  of  the  variance  each  time  the  aircraft 
passes  on  this  side.   A  second  reason  for  a  change  in 
variance  is  the  aircraft's  location  in  the  pattern. 
However,  these  modulations  are  not  as  apparent  since  the 
range  variation  dominates.   RMS  errors  were  obtained  by 
repeating  these  runs  20  times  and  are  shown  in  Figures  16 
thru  18.   The  fluctuations  generated  by  the  flight  path  are 
apparent,  more  so  in  the  six-state  system  than  in  the 
two-state  system.   Also,  steady  state  values  cf  the  error 
gradually  increase.   This  is  true  because  the  average  range 
to  the  so  no  buoy  increases  as  the  buoy  moves  farther  away 
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from  the  center  of  the  pattern.   (It  will  be  demonstrated 
later  that  the  errors  are  range  dependent.)   The  covariance 
of  the  six-state  system  increases  from  400  yards  at  25 
minutes  to  800  yards  at  the  end  of  the  simulation.   The 
two-state  system  increases  slowly  from  250  yards  until  it  is 
hammered  by  the  Schuler  cycle.   It  is  worthy  of  note  that 
vithout  the  navigational  errors  stemming  from  the  Schuler 
cycle  the  RMS  errors  of  the  two-state  system  are 
considerably  lover  than  those  of  the  six-state  system. 

Table  1  compares  the  RMS  errors  of  these  two  systems  to 
the  results  obtained  by  Orincon  for  their  13-state  system. 
(The  values  for  the  two-state  and  six-state  systems  are 
taken  from  the  covariance  at  a  point  35  minutes  into  the 
simulation.)   The  first  line  considers  measurements  of 
bearing,  doppler  velocity,  and  drift  angle.   Navigational 
errors  are  present  including  a  Schuler  cycle.   As  expected 
there  is  an  increase  in  relative  position  error  because  of 
the  reduction  from  13  to  6  states.   The  next  line  considers 
the  effect  of  the  Jlooal  Positioning  System  (3-PS),  with  very 
accurate  aircraft  positioning  information,  on  the  relative 


13-state 

6-sta te 

2-s tate 

Schuler 
SRS-dv*DA 

394 

453 

No  Schuler 
SHS+DV+DA+GPS 

450 

200 

Table  1.   Comparison  cf  RMS  results  in  yards. 
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position  errors  of  the  sonotuoy.   In  the  13-state  system 
Orincon  observed  that  the  errors  in  the  states  pertaining  to 
the  aircraft  inertials  improved  markedly;  relative  position 
errors  most  likely  remained,  the  same.   This  is  exactly  vhat 
was  found  in  the  analysis  of  the  six-state  system:  that  the 
relative  position  errors  were  unaffected  by  navigational 
errors.   However,  the  two-state  system  is  very  much  affected 
by  navigational  errors.   Relative  position  errors  were 
significantly  less  for  the  two-state  system  when  accurate 
aircraft  positions  were  known. 

To  obtain  a  more  realistic  analysis  of  their  performance, 
each  algorithm  was  tested  using  the  circular  and  square 
flight  patterns  of  Figures  6  and  7.   The  results  of  the 
first  simulation  run  for  the  circular  path  at  a  range  cf  15 
NM  are  shown  in  Figures  19  thru  21.   In  this  case,  the  v.avy 
nature  of  the  covariance  is  due  to  the  location  of  the 
aircraft  in  the  pattern  since  the  range  to  the  sonotuoy 
remains  relatively  constant.   For  instance,  as  the  aircraft 
passes  directly  north  or  south  of  the  sonobuoy  the  ability 
to  correct  E-W  errors  is  greatest.   Therefore,  the  E-W 
variance  reaches  a  minimum  value  at  this  time.   The  N-S 
variance  is  affected  in  a  similar  manner.   The  two-state 
filter  is  not  as  sensitive  to  this  as  the  six-state  filter. 
Also,  the  effect  of  the  Schuler  cycle  on  the  two-state 
filter  can  be  seen  in  Figure  20. 

The  RMS  errors  are  shown  in  Figures  22  thru  24.   The 
steady  state  RMS  value  of  error  predicted  by  the  filter  is 
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425  yards  for  the  two-state  and  750  yards  for  the  six-state 

system.   The  actual  RMS  error  in  the  relative  position  of 

the  sonobuoy  is  250  to  1200  yards  for  the  two-state  and  400 

yards  for  the  six-state  system.   Without  Schuier  cycle  error 

the  two-state  system  drops  to  a  steady  250  yards,  better 

than  the  six-state  system.   It  is  interesting  to  note  that 

the  measured  errors  are  significantly  less  than  those 

predicted  ty  the  filter.   Both  measured  values  of  RMS  error, 

a  and  a   ,  are  in  close  agreement;  that  is,  the  deviation  of 
x      a 

the  filter's  estimated  position  about  its  mean  is  generally 
the  same  as  the  deviation  of  the  estimated  position  about 
the  actual  relative  position  of  the  sonobuoy.   The  closer 
these  values  are  to  one  another,  the  more  confidence  can  te 
placed  in  the  analysis.   The  only  exception  to  this  is 
Figure  23  which  shows  the  reaction  of  the  two-state  filter 
to  the  Schuier  cycle.   One  complete  cycle  with  a  period  of 
84  minutes  is  obvious.   The  Kalman  filter  does  not  recognize 
this  cycle  since  the  modeling  equations  do  not  account  for 
it.   The  dip  in  the  center  of  the  two  peaks  is  once  again 
the  result  of  flight  path  geometry.   The  aircraft  is  in  such 
a  position  relative  to  the  sonobuoy  that  the  measurements 
provide  enough  information  to  correct  for  Schuier  cycle 
errors.   However,  it  is  not  in  this  position  long  enough  to 
influence  the  errors  anymore  than  it  does. 

The  same  analysis  was  performed  using  the  square 
pattern  of  Figure  7.   This  scenerio  has  mere  drift  than  the 
other  and  in  a  different  direction  providing  the  algorithms 
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with  another  motion  to  track.   The  outcome  was  generally  the 
same  as  can  be  seen  in  Figures  25  thru  30.   There  is  no 
measurable  difference  in  the  results  from  the  tvo  different 
flight  paths.   Flying  a  straight  path  does  not  adversely 
affect  the  results  as  was  concluded  by  the  Oricon 
Corporation  for  their  13-state  system.   The  wavy  nature  of 
the  covariance  is  again  a  result  of  the  flight  path.   It  is 
not  as  smooth  as  before  because  a  square  pattern  was  flown 
instead  of  a  circular  one.   Another  look  at  errors  caused  in 
the  two-state  system  by  the  Schuler  cycle  can  be  seen  in 
Figure  30.   The  algorithm  does  make  some  corrections  to 
these  errors  but  they  are  not  as  effective  as  before.   In 
this  case  the  time  the  aircraft  was  in  a  position  to  ma^ce 
the  corrections  aid  not  coincide  v,ith  the  time  the  peak 
errors  occurred.   There  appears  to  be  no  way  to  predict  the 
optimum  time  and  place  for  the  aircaft  to  te  without  knov.ing 
when  and  hew  the  Schuler  cycle  is  occurring. 

The  distance  the  aircraft  is  from  the  sonotuoy  is 
directly  related  to  the  accuracy  the  algorithms  can  achieve. 
Simulations  were  performed  at  ranges  of  5,  15,  33,  and  45 
nautical  miles  using  a  circular  flight  path  with  the 
sonobuoy  in  the  center.   (The  RMS  plots  can  be  fcund  in 
Appendix  G.)   The  steady  state  errors  were  observed  and  are 
plotted  in  Figure  31.   Solid  lines  represent  the  RMS  values 
of  the  covariance  and  dashed  lines  represent  the  actual 
errors.   In  all  cases  the  errors  increase  with  distance 
>.hich  is  intuitively  satisfying.   Mathematically,  the 
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covariance  increases  because  the  H  matrix  becomes  smaller 
with  an  increase  in  range.   Consequently,  the  updated  value 
of  the  covariance,  P(+),  is  larger  from  Equation  (19).   For 
the  six-state  system  there  is  an  increase  in  actual  error 
from  240  yards  at  5  NM  to  750  yards  at  45  NM .   And  for  the 
two-state  system,  actual  errors  increase  from  13*3  yards  at  5 
NM  to  500  yards  at  45  NM.   (Note  that  for  the  Uo-state 
system,  Schuler  cycle  errors,  indicated  by  circles,  peak 
approximately  1300  yards  above  these  plotted  values!  the 
lower  values,  indicated  by  x's,  were  used  since  they  were 
more  reliable  for  comparison  between  ranges.)   There  is  a 
slight  decrease  in  the  slopes  of  all  the  error  curves  as 
distance  increases.   Hovever,  this  decrease  is  small  and  the 
curves  might  well  have  been  interpeted  as  linear  within  the 
limits  of  the  analysis. 

The  frequency  at  which  measurements  are  made  also 
affects  the  accuracy  of  the  estimate.   Each  algorithm  \>as 
tested  with  measurement  intervals  of  4,  10,  20,  and  30 
seconds  on  one  sonobuoy.   The  outcomes  of  the  steady  state 
RMS  errors  are  plotted  in  Figure  32.   (Once  again  the 
two-state  system's  errors  do  not  show  the  effect  of  the 
Schuler  cycle.)   The  actual  errors  in  the  six-state  system 
increase  steadily  as  the  measurement  interval  increases. 
However,  the  covariance  decreases  rapidly  at  first  to  a 
minimum  value  somewhere  around  13  seconds  and  then  increases 
with  increasing  interval.   An  abbreviated  run  with  a  two 
second  measurement  interval  confirmed  these  results. 
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Propagation  effect 


Observability  effect 


Figure  33.   Factors  which  influence  the  covariance. 

Figure  33  helps  to  explain  this  outcome.   There  are  tv.o 
conditions  which  affect  the  covariance.   The  most  obvious  is 
the  increase  in  covariance  due  to  an  increase  in  the 
propagation  interval.   This  occurs  because  the  plant  matrix, 
vhich  is  a  function  of  time,  affects  the  covariance  in 
Equation  (16).   The  second  condition  is  a  consequence  of  the 
model's  observability.   If  the  relative  bearing  betv-een  the 
aircraft  and  the  scnobuoy  does  not  change  then  the  state  is 
not  observable.   In  other  words,  two  successive  bearings 
must  intersect  to  determine  an  estimated  position.   By 
allowing  the  interval  between  measurements  to  become  too 
small,  the  aircraft  is  unable  to  make  a  significant  change 
relative  to  the  sonobuoy.   The  covariance  increases  as  the 
conditions  approach  those  which  maKe  the  problem 
unobservable .   There  appears  to  be  an  optimum  frequency  with 
which  to  make  measurements  en  one  buoy,  namely  10  seconds 
for  this  range  of  15  NM  and  speed  of  180  knots. 


6^+ 


The  two-state  system  exhibited  similar  characteristics 
in  Figure  32.   In  this  case  the  covariance  never  does 
increase  in  the  range  of  intervals  studied.   The  actual 
errors  remain  constant.   Since  the  plant  matrix  fcr  this 
system  is  equal  to  the  identity  matrix  and  the  Q  matrix  is 
prescribed  to  be  constant,  the  propagation  interval,  &t, 
does  not  effect  the  covariance.   However,  the  effects  of  the 
observability  conditions  do  cause  the  covariance  to  increase 
vhen  At  becomes  too  small.   It  is  believed  that  as  At 
becomes  smaller  the  actual  errors  would  begin  to  increase 
also  . 

The  statistical  plots  for  the  range  tests  found  in 

Appendix  B  show  that  the  Schuler  cycle  causes  larger  peaks 
in  the  actual  error  as  the  range  from  the  sonobuoy 
increases.   However,  the  peaks  are  decreasing  with  time. 
Further  simulation  revealed  that  the  peak  errors  decrease  to 
a  steady  value  about  10J0  yards  greater  than  non-Schuler 
cycle  errors.   This  coincides  with  the  amplitude  of  the 
Schuler  cycle  as  it  was  programmed  for  this  particular 
flight  path.   Therefore,  in  the  steady  state  a  Schuler  cycle 
may  cause  additional  estimation  errors  approximately  equal 
to  its  amplitude. 

Altitude  was  also  tested  for  any  effect  it  might  have  en 
these  algorithms.   Using  the  circular  pattern  at  a  range  of 
15  MM,  altitudes  of  300,  3000,  10,000,  and  20,000  feet  were 
tested.   Neither  the  six-state  nor  the  two-state  system 
shoved  any  effect  for  these  changes  in  altitude.   It  is 


65 


believed  that  at  very  small  ranges  (less  than  two  nautical 
miles)  altitude  might  hamper  good  position  estimation. 
However,  for  the  majority  of  the  time  altitude  is  of  no 
concern . 
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Figure  13.   Est.  position  and  relative  errors  for  the  six-state  system  using 
Orincon's  pattern.   DR=(a/c  225-5.2  +  Schuler,  b  045-5.5),  At=20s,  Alt=3000'. 
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Figure  14.   Est.  position  and  relative  errors  for  the  two-state  system  using 
Orincon's  pattern.   DR=(a/c  225-5.2  +  Schuler,  b  045-5.5),  At=20s,  Alt=3000'. 
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Figure   15.      Est.    position  and  relative   errors   for   the   two-state   system  using 
Orincon's   pattern.      DR=(a/c   225-5.2  +  No   Sen.,   b   045-5.5),    At=20s,   Alt=3000'. 
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Figure  19.   Est.  position  and  relative  errors  for  the  six-state  system  using 
circular  pattern  at  15NM.   DR=(a/c  180-1  +  Schuler,  b  270-1),  At=20s,  Alt=3000' 
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Figure  20.   Est.,  position  and  relative  errors  for  the  two-state  system  using 
circular  pattern  at  15NM.   DR=(a/c  180-1  +  Schuler,  b  270-1),  At=20s,  Alt=3000' 
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Figure  21.   Est.  position  and  relative  errors  for  the  two-state  system  using 
circular  pattern  at  15NM.   DR=(a/c  180-1  +  No  Sen.,  b  270-1),  At=20s,  Alt=3000' 
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Figure  25.   Est.  position  and  relative  errors  for  the  six-state  system  using 
square  pattern  at  15NM.   DR=(a/c  000-2.5  +  Schuler,  b  090-2),  At=20s,  Alt=3000' 
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Figure  26.   Est.  position  and  relative  errors  for  the  two-state  system  using 
square  pattern  at  15NM.   DR=(a/c  000-2.5  +  Schuler,  b  090-2),  At=20s,  Alt=3000' 
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Figure  27.   Est.  position  and  relative  errors  for  the  two-state  system  using 
square  pattern  at  15NM.   DR=(a/c  000-2.5  +  No  Sch. ,  b  090-2),  At=20s,  Alt=3000 ' 
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V.   CONCULSION 

The  following  are  the  conculsions  reached  at  the  end  of 
this  research: 

1.  The  six-state  algorithm  developed  in  this  report 
produced  estimates  of  relative  sonobuoy  position  with  a 
standard  deviation  in  the  error  of  500  yards  at  a  range 
of  20  NM.   Aircraft  drift,  which  included  the  Schuler 
cycle,  had  no  effect  on  this  system  since  the  states 
were  defined  relative  to  the  aircraft. 

2.  The  two-state  algorithm  developed  in  this  report 
produced  estimates  of  relative  sonobuoy  position  with  a 
standard  deviation  of  300  yards  plus  an  amount  equal  to 
the  amplitude  of  the  Schuler  cycle  error  at  a  range  of 
20  NM .   When  nc  Schuler  cycle  navigational  errors  were 
present,  the  two-state  system  cut-performed  the 
six-state  system. 

3.  Both  algorithms  shoved  an  increase  in  steady  state  RMS 
error  with  an  increase  in  range. 

4.  The  six-state  system  showed  an  increase  in  actual 
steady  state  RMS  error  for  an  increase  in  the  time 
interval  between  measurements?  however,  the  covariance 
of  this  system  indicated  minimum  error  occurred  with  a 
10  second  interval  on  one  sonobuoy  at  a  range  of  15  NM 
and  an  aircraft  speed  of  180  lets.   (This  translates  to  a 
lh   degree  change  in  relative  bearing.)   The  two-state 
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system  shewed  no  change  in  the  actual  steady  state  RMS 
error  for  intervals  between  4  and  33  seconds  and  only  a 
slight  increase  in  these  errors  as  a  60  second  interval 
was  approached.   The  covariance  of  the  system  actually 
decreased  with  increasing  interval. 

5.  Neither  system  showed  any  dependence  on  altitudes 
within  the  operational  limits  at  a  range  of  15  MM. 

6.  No  adverse  effect  was  observed  while  flying  straight 
tracks  as  opposed  to  curved  ores. 

7.  Steady  state  RMS  error  was  affected  by  two  conditions 
as  measurement  intervals  changed.   The  longer  the 
interval  the  more  the  error  tended  to  increase, 
especially  for  the  six-state  algorithm,  as  a  result  of 
system  propagation.   The  shorter  the  interval  the  more 
the  error  tended  to  increase  as  the  systems  encroached 
upon  their  non-observability  condition.   An  optimum 
measurement  interval  existed  for  each  system  which  was  a 
function  of  range. 
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VI-   SUMMARY 

Both  the  six-state  and  the  tv,o-state  algorithms  were 
atle  to  successfully  track  the  sonobuoy.   However,  the 
six-state  system  was  the  best  choice  given  the  navigational 
system  in  use  by  the  Orion.   It  was  capable  of  providing 
estimated  sonobuoy  positions  at  a  range  of  22  NM  such  that 
the  errors  had  a  standard  deviation  of  502  yards.   The 
computer  time  required  to  process  one  measurement  was  on  the 
order  of  222  milliseconds.   On  the  other  hand,  the  much 
simpler  two-state  system  required  only  2S  milliseconds  or 
one  seventh  (y)  the  amount  of  computer  time.   However,  its 
estimate  of  sonobucy  position  was  degraded  by  Schuler  cycle 
errors  in  the  navigational  system.   Without  these,  this 
system  was  subject  to  errors  having  a  standard  deviation  of 
322  yards  at  a  range  of  20  NM.   The  Schuler  cycle 
periodically  increased  these  errors  by  an  amount  equal  tc 
its  amplitude. 

The  algorithms  i*ere  very  much  affected  by  the  aircraft's 
flight  path  including  the  relative  location  of  the  sonobucy, 
the  range  to  the  sonobuoy,  and  the  amount  of  angular  change 
in  the  bearing  to  the  sonobuoy  as  measurements  v,ere  taken. 
Since  the  sonobuoy  position  was  to  be  estimated  in 
tv»  o-dimensions ,  the  aircraft  was  required  tc  maneuver  in 
such  a  vay  so  as  to  provide  information  in  both  dimensions. 
For  example,  as  the  aircraft  crossed  north  or  south  cf  the 
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sonobuoy  the  E-W  errors  decreased  while  the  N-S  errors 
increased.   The  trend  reversed  when  the  aircraft  flew  east 
or  vest  of  the  sonobuoy.   It  \»as  also  observed  that  steady 
state  RMS  errors  grew  with  range  at  a  slightly  decreasing 
rate.   They  were  twice  as  large  at  45  NM  as  they  were  at  15 
NM.   Finally,  the  sonobuoy  position  was  not  observable  if 
the  relative  bearing  to  the  sonobuoy  did  not  change  implying 
that  two  lines  of  bearing  must  cross  in  order  to  obtain  an 
estimated  position.   If  the  amount  of  angular  change  between 
each  measurement  was  too  small  it  began  to  affect  the  errors 
because  this  condition  was  being  approached.   Anytime  the 
aircraft  flew  in  such  a  way  as  to  reduce  the  rate  of  bearing 
change,  or  the  frequency  at  which  measurements  were  made  >as 
too  high,  this  condition  was  prevalent.   It  was  observed 
that  at  least  1%  degrees  of  change  was  required  to  avoid 
this  problem. 

The  sii-state  system  showed  no  tendency  toward 
divergence  as  long  as  initial  conditions  were  reasonable. 
The  two-state  system,  on  the  other  hand,  did  have  a  tendency 
to  diverge  when  the  gains  on  the  drift  filter  were  not 
adjusted  correctly.   For  example,  if  the  gain  was  to  high 
the  estimated  drift  was  heavily  dependent  on  the  last  fev 
measurements  and  was  too  quickly  affected  by  changes  in  the 
position  filter.   But,  it  was  also  necessary  for  the  drift 
filter  to  have  non-zero  gains  in  the  steady  state  in  order 
to  provide  some  insurance  against  an  initially  wrong 
estimate  of  sonobuoy  drift.   With  these  modifications  maae 
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the  algorithm  operated  correctly  during  all  tests  and,  in 
fact,  had  less  error  in  some  cases  than  did  the  six-state 
system. 

It  is  the  opinion  of  the  author  that  either  one  of  these 
algorithms  is  an  improvement  on  the  nistorical  method  of 
"mariting-on-top"  the  sonobuoy.   further,  it  adds  a  great 
deal  mere  uncertainty  as  to  exactly  ho*  much  accuracy  could 
be  obtained  in  the  old  way.   Non-linear  navigational  effects 
were  not  even  considered  and  the  number  of  measurements  made 
to  determine  sonobuoy  position  was  on  the  order  of  10  or 
less  for  any  one  sonobuoy  pattern.   Nevertheless,  in  an 
operational  sense  the  accuracies  achieved  here  would  need  to 
be  improved.   Without  the  non-linear  navigational  errors  the 
accuracies  could  have  been  much  better.   Eowever,  the 
results  of  this  research  do  provide  another  perspective  for 
the  problem  at  hand. 
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APPENDIX  A 
THE  COMPUTER  PROGRAMS 


MAIN 


Simulation  setup 


Do  for  the  no. 
of  simulations 


Read  measurement i 
cos  &    and  antenna 
doppler  velocity  and  drift  angle-* 
bouy,  time,  a/c  lat .  and  long, 
altitude,  heading,  pitch,  roll 
accelerometer  outputs* 


YES 


Initialize  buoy's 
state  vector  and 
covariance  matrix 


YES 


Process  results 


Retrieve  buoy's- 
state  vector 
covariance  matrix 
time  of  measurement 


(Stop) 


Convert  lat/long  to  NM 
Determine  antenna  baseline  vector 
At  =  present  time  -last  time 
Z3  =  heading  -  drift  angle* 


CALL 
FILTER 


Store  new- 
state  vector 
covariance  matrix 
time  of  measurement 


Process  simulation  results 


Figure  Jk- .      Flow  chart: 
main  program 


*  Required   for  six-state 
system   only. 
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SUBROUTINE  FILTER 
SIX -3 TATE 


Calculate 
plant  matrix  ($. 


constant  matr 


ix  A,,. 


Propagate 
state  vector 

\  -  h-A-i  +  K 

covariance  matrix 


3k  aK*ic-i<-isk-AVVl]* 


Do  1=1.3 


YES 


z.    -   cos  a 
Calculate 
dz  =  z,  -  z. 
H  =  Hj 
R  =  R, 


no 


YES 


'2     d 
Calculate 

dz  =  z~  - 

H  =  H2 

R  =  R, 


z3  B  He 
Calculate 

dz  =  z,  - 

H  =  H.. 

R  =  R, 


Update 
gain 

-T  T 

f  -  vv 

state  vector 


A   -    f  f 


t<   k 


k    k 

covariance 

\  ■  sk  ci 

T 
P      =    SS 


+   G   dz 


T    J 
ff    ]* 


MO 


YES 


Xk    "    Xk 


Sk    " 


_(  Next  ) 


© 


Figure  35 


Flow  chart 
six-state 


Subroutine  FILTER 
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SUBROUTINE  FILTER 
TWO -STATE 


Save  last  position  & 

velocity  +     . 

Xk-f  Uk-1 

Propagate 

state  vector 

Xk  -   Xk-1  +  °k-lAt 

covariance 

sk =  K-isk-i +  Qk-i]i 

Calculate 

X  "  Xk  +  Xa 

H(l,l)  =  Htl 

H(l  ,2 

)  =  H12 

Update 

gain 

-T   T           T                t 
f  =  S^H^     C<  =    fXf  +  R     G  =  ~ 

V 

state  vector 

dz  =  a  -  z           z  =  cose 

x;  =  x-  ♦  g  dz 

covariance 

<=3ktl-i  ffT] 

P   =  SST 

Drift 

measurement 

Zu  =  (x£  -  X^jJ/at  ,   limit  to  20 

kts 

gain 

C  = 

i 

limi 

1-  tn     ^     r 

linimum 

n  +  1 

34  min 

inhibit  first  5  dr 

ift  calculations 

(C  =  0) 

n  = 

1+1 

Drift 

U*  =  U*   +  C  (Z   -  V?  , 
k    k-1       u    k-1 

0 


Figure  36. 


Flow  chart:      Subroutine 
FILTER    two-state 


9/4. 


GLOSSARY  OF  COMPUTER  VARIABLES 
SIMULATION  PROGEAM 

A       Vector  of  aircraft  accelerations 

ALPHAI  Result  of  -^— ^ or  - 

fTf  +  R 
ALT     Altitude 

ANGLE    Direction  of  sonobuoy  drift  from  North 

ANT     Antenna 

AVEL     Initial  aircraft  velocity 

T 
AS      Matrix  result  of   I  -  ^~ 


B  Matrix  result  of  (eSS1  +  Q> 

BOUY  Sonobuoy  RE  channel  or  number 

C  Gain  of  drift  filter 

CA  Result  Of  At/(84  min.) 

C^  Constant  matrix  in  six-state  system's  equations 

CKL  Rate  of  change  of  nautical  miles  for  a  change  in 
degrees  of  longitude 

D  Array  which  holds  all  items  requiring  storage  from 
one  measurement  to  the  next 

CELT  Interval  between  measurements,  At 

DELX1  Relative  position,  axx 

DELX2  Relative  position,  ax2 

DELX3  Relative  position,  ax  3 

DRIFT  Magnitude  of  sonobuoy  drift  in  yards 

DSEED  Seed  for  the  random  number  venerator 

DTIMS  Array  of  initial  drop  times  for  the  sonobuoy 

DZ  Res  idual ,  dz 

SDELX1  Error  in  the  est.  relative  position  of  sonobuoy,  &Lxl 
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ELELX2  Error  in  the  est.  relative  position  of  sonobuoy,  ax 2 

EOSSO  Matrix  result  of  e  $TSTS  23 

EPSLN  Age  weighting  factor, 

ESST  Matrix  result  of  e  S  5T 

EX1  Error  in  est.  geographical  position  of  sonobuoy, a*! 

EX2  Error  in  est.  geographical  position  of  sonobuoy, Ax2 

F  f  vector  from  Carlson  square  root  technique 

T 

FFT     Matrix  result  of  ff 

T 
ff 


EFTA  Matrix  result  of 

FT  Transpose  of  f 

3-  Gain  matrix 

GAMMA  Aircraft  pitch  angle,  y 

GDZ  Vector  result  of  G  dz 

GR  Gamma  in  radians 

H  Aircraft  heading 

H  The  H  matrix  when  used  in  subroutine  FILTER 

HOUR  Time  in  hours 

HR  Heading  in  radians 

HT  Transpose  of  the  H  matrix 

I  Index 

I  Identity  matrix  when  used  in  subroutine  FILTER 

IEND  Indicator  of  EOF 

T 
ff 

IFFTA  Matrix  result  of   I  -  —^- 

IN  File  from  which  data  is  read 

INRES  File  from  which  true  scnobuoy  positions  are  read 

IOUT  Device  on  which  print  is  made 

IOUTF  File  in  which  results  are  written 
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IOUTZ  File  in  which  random  vectors  are  inspected 

IRUN  Index  for  run  number 

IZ  Index  for  measurement,  z  =  1,2,3 

K  Incremental  measure  of  time 

LAT  Latitude 

LATI  The  latitude  of  the  origin  of  the  local  earth  fixed 
coordinate  system 

LONG  Longitude 

LONSI  The  longitude  of  the  origin  of  the  local  earth  fixed 
coordinate  system 

MARK  Indicates  whether  information  is  to  be  stored 
or  retrieved 

MEAM  Mean  estimated  sonobuoy  position,  mi 

MEAN2  Mean  estimated  sonobuoy  position,  u 2 

MIN  Time  in  minutes 

N  Total  number  of  measurements  made  on  all  sonobuoys 

N  Index  *hen  used  in  subroutines  IN  ITAL  and  DATA 

MRUNS  Integer  equivalent  of  RUNS 

0  Plant  matrix,  $ 

OS  Matrix  result  of  $  S 

OST  Transpose  of  $  ST 

OSSO  Matrix  result  of  $TSTS  $ 

OX  Matrix  result  of  §   X 

?  Covariance  matrix 

PHI  Aircraft  roll  angle,  f 

PI  TT 

PR  Phi  in  radians 

Q  Variance  cf  plant  noise  from  N(0fQ) 

R  Variance  of  measurement  ncise  from  N(3,R) 
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RAD 
RANGE 
RANSE2 
RBI 

RB2 
RB3 
RBX 
RNO 

RUNS 
S 

SI 

S2 

SEC 

SE 

SIGMA 

SIGMA? 

SIGMAX 

SKL 

SNSXT 

SST 

ST 

SUMP1 

SUMP2 

SUMSQ1 

SUMSQ2 

SUMX1 


tt  /180 

Range  to  sonobuoy  from  aircraft 

Squared  value  of  range 

Unit  vector  component  of  antenna  baseline 

Unit  vector  component  of  antenna  baseline 

Unit  vector  component  of  antenna  baseline 

Result  of  (RB1  &xx     +   RB2ax2  +  RB3ax3) 

Normally  distributed  random  number  from  N ( d  ,  1 ) 

Number  of  simulation  runs  to  be  made 

Array  in  Carlson's  square  root  technique  which 
represents  square  root  of  covariance  matrix 


Square  root  of  variance,  /Pn 

Square  root  of  variance,  ^fiz 

Time  in  seconds 

Matrix  result  of  S  f 

RMS  error  in  the  estimate  about  the  actual  location 


RMS  value  from  the  Kalman  filter,  /Pj  \  +  F22 

RMS  error  in  the  estimate  about  the  mean  estimate 

Rate  of  change  of  nautical  miles  for  a  change  in 
degrees  of  latitude 

Propagated  matrix  of  S 

Matrix  result  of  S  ST 

Transpose  of  S  matrix 

Array  for  the  sum  of  P(l,l 

Array  for  the  sum  of  P(2,2 


over  NRUNS 
over  NRUNS 


Array  for  sum  of  ax  over  NRUNS 
Array  for  sum  of  ax  over  NRUNS 
Array  for  sum  of  a'x.  over  NRUNS 
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SUMX2 


Array  for  sum  of  Ax2  over  NRUNS 


TIME 

0 

Ul 

U2 

VAfil 

VAR2 

VARX1 

VARX2 

VELCT 
VELCT2 

x 

XAI 

XAT 

XI 

X2 

X3 

XB1 

X32 

XB1I 

XB2I 

XNEXT 

Z 


Last  time  in  seconds  the  sonobuoy  was  processed  Dy 
the  filter 

Time  from  aircraft  clock  in  seconds 

Sonobuoy  drift  vector 

Old  value  of  U(l) 

Old  value  of  U(2) 

Variance  of  the  estimated  sonobuoy  position  about 
the  actual  location   af 

Variance  of  the  estimated  sonobuoy  position  about 
the  actual  location  o\ 

Variance  of  the  estimated  sonobuoy  position  about 

the  mean  estimated  oosition,  a2, 

xl 

Variance  of  the  estimated  sonobuoy  position  about 
the  mean  estimated  position,  a2 

Relative  velocity,  aV 

Squared  value  of  relative  veclocity 

State  vector 

Aircraft  inertial  position,  X 

a 

Aircraft  true  position,  X 

a 

Old  value  of  X(l) 
Old  value  of  X(2) 

Vector  including  XB1,  XB2,  velocity  of  sonotuoy,  and 
two  zeros  used  in  the  six-state  system 


True  geographical  location  of  sonobuoy,  x 
True  geographical  location  of  sonobuoy,  x 
Array  of  initial  sonobuoy  positions 
Array  of  initial  sonobuoy  positions 
Propagated  value  of  X 
Measurement  vector 
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ZHAT    Estimate  of  z 

ZN      Normally  distributed  random  number  from  N(z,R) 

ZU      Vector  of  velocity  measurement  for  drift  filter 
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APPENDIX  B 
COMPUTER  PROGRAM  FOR  DATA  GENERATION 


DATA  GENERATION 


Set  initial  conditions: 

sonobouy  position  and  velocity 
aircraft  position  and  velocity 
aircraft  heading 
delt 


LINE 


Set: 

Phi  ■  0 
Accel  -  0 

Calculate: 
aircraft  position 

Call 
FILE 


0 


CURVE 


Calculate : 

Phi 
Accel 


Calculate : 
aircraft  position 


FILE 


Call 

TI 

ME 

Calcul 

ate : 

aircraft  drift 

sonobouy  number 

sonobouy  position 

measu 

rements,  I 

Set: 

gamn 

la  «  0 

Write 

data  into 

storage  file 

Call 
FILE 


© 


0 


Figure   37.      Flow  chart:      data   generation  program. 
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GLOSSARY  01  COMPUTER  VARIABLES 
DATA  GENERATION  PROGRAM 

A  Vector  of  aircraft  accelerations 

AL  Altitude  in  nautical  miles 

ALT  Altitude 

ANT  Antenna 

AVEL  Initial  aircraft  velocity 

AX1  Amplitude  of  Schuler  cycle,  A,. 

AX2  Amplitude  of  Schuler  cycle,  A  ^ 

3  Intercept  on  N-S  axis  of  straight  track 

BANG  Sonotuoy  drift  direction 

BOUY  Sonobucy  RE  channel  cr  number 

BVEL  Sonobuoy  drift  direction 

CKL  Rate  of  change  of  nautical  miles  for  a  change  in 
degrees  of  longitude 

D  Length  of  straight  track 

DELT  Interval  between  measurements,  At 

DELV1  Relative  velocity,  Avj 

DELV2  Relative  velocity,  Av 2 

DELX1  Relative  position,  ^xl 

DELX2  Relative  position,  ax 2 

DH  Total  change  in  heading  while  on  curved  track 

DSEEI)  Seed  for  the  random  number  generator 

END  Indicator  of  EOE 

i  Gravitational  acceleration 

GAMMA  Aircraft  pitch  angle,  y 
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GR  Gamma  in  radians 

H  Aircraft  heading 

KNEW  Final  heading  on  curved  track 

HOUR  Time  in  hours 

HR  Heading  in  radians 

IOUT  Device  on  \*hich  print  is  made 

IOUTP  File  in  which  results  are  written 

IOUTR  File  in  which  true  sonobuoy  positions  are  written 

IOUTZ  File  in  fchich  random  vectors  are  inspected 

ITIME  Time  from  aircraft  clock  in  seconds 

KXl  Rate  of  linear  navigational  drift,  K 

KX2  Rate  of  linear  navigational  drift,  K 

LAT  Latitude 

LATI  The  latitude  of  the  origin  of  the  local  earth  fixed 
coordinate  system 

LONG  Longitude 

LONGI  The  longitude  of  the  origin  of  the  local  earth  fixed 
coordinate  system 

LTIME  Array  of  initial  drop  times  for  the  sonobuoy 

M  Slope  of  straight  track 

MARK  Indicates  the  first  time  a  buoy  is  to  be  proccessed 

MIN  Time  in  minutes 

N  Index 

NBOUY  Total  number  of  sonobuoys 

NSTOP  End  of  line  or  curve  series 

PHI  Aircraft  roll  angle,  ft 

PI  TT 

FI2  2^ 
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PR  Phi  in  radians 

R  Radius  of  curved  track 

RBI  Unit  vector  component  of  antenna  baseline 

RB2  Unit  vector  component  of  antenna  baseline 

RB3  Unit  vector  component  of  antenna  baseline 

S  Distance  aircraft  travels  in  time  At 

SEC  Time  in  seconds 

SKL     Rate  of  change  of  nautical  miles  for  a  change  in 
degrees  of  latitude 

TRACK  Direction  of  aircraft  travel  over  the  ground 

U  Direction  to  center  cf  curved  track 

V  Direction  from  center  of  curved  track  to  aircraft 

VA1  Velocity  of  aircraft 

VA2  Velocity  of  aircraft 

V31  Sonobuoy  drift 

VB2  Sonobucy  drift 

WANG  Wind  direction 

WVEL  Wind  velocity 

XAI     Aircraft  inertial  position,  X 

a 

XAT     Aircraft  true  position,  X 

i 

i     nf    <;oTinhiinv.     x 

1 


a 
XB1     True  geographical  location  of  sonobuoy,  x 


XB2  True  geographical  location  of  sonobuoy,  x 

XC1  Location  of  the  center  of  the  curved  track 

XC2  Location  of  the  center  of  the  curved  track 

XXI  Incremental  change  in  X 

XX2  Incremental  change  in  X 

2 

Z  Measurement  vector 
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APPENDIX  C 
Range  Analysis 
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Figure  38.   RMS  errors  for  two-state  system  (top)  and  six- 
state  system  (bottom)  using  circular  pattern  at  5  NM. 
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Figure  39.   RMS  errors  for  two-state  system  (top)  and  six- 
state  system  (bottom)  using  circular  pattern  at  30  NM. 
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Figure  ^0 .   RMS  errors  for  two-state  system  (top)  and  six- 
state  system  (bottom)  using  circular  pattern  at  k$   NM. 
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Figure  4-1.   RMS  errors  for  two-state  system  (top)  and  six- 
state  system  (bottom)  using  square  pattern  with   t=4  sec. 
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Figure  42.   RMS  errors  for  two-state  system  (top)  and  six- 
state  system  (bottom)  using  square  pattern  with   t=10  sec 
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Figure  kj.      RMS  errors  for  two-state  system  (top)  and.  six- 
state  system  (bottom)  using  square  pattern  with  t=30  sec 
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Altitude  Analysis 
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Figure  ^3  •   RMS  errors  for  two-state  system  (top)  and  six- 
state  system  (bottom)  using  circular  pattern,  Alt=300'. 
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Figure  4-5.   RMS  errors  for  two-state  system  (top)  and  six- 
state  system  (bottom)  using  circular  pattern,  Alt=10,000'. 
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Figure  46.   RMS  errors  for  two-state  system  (top)  and  six- 
state  system  (bottom)  using  circular  pattern,  Alt=20,000' 
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